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ABSTRACT 

When our Sun was young it rotated much more rapidly than now. Observations of young, rapidly rotating 
stars indicate that many possess substantial magnetic activity and strong axisymmetric magnetic fields. We 
conduct simulations of dynamo action in rapidly rotating suns with the 3-D MHD anelastic spherical harmonic 
(ASH) code to explore the complex coupling between rotation, convection and magnetism. Here we study 
dynamo action realized in the bulk of the convection zone for two systems, one rotating at three times the 
current solar rotation rate and another at five times the solar rate. We find that substantial organized global- 
scale magnetic fields are achieved by dynamo action in these systems. Striking wreathes of magnetism are 
built in the midst of the convection zone, coexisting with the turbulent convection. This is a great surprise, 
for many solar dynamo theories have suggested that a tachocline of penetration and shear at the base of the 
convection zone is a crucial ingredient for organized dynamo action, whereas these simulations do not include 
such tachoclines. Some dynamos achieved in these rapidly rotating states build persistent global-scale fields 
which maintain amplitude and polarity for thousands of days. In the case at five times the solar rate, the dynamo 
can undergo cycles of activity, with fields varying in strength and even changing polarity. As the magnetic fields 
wax and wane in strength, the primary response in the convective flows involves the axisymmetric differential 
rotation, which begins to vary on similar time scales. Bands of relatively fast and slow fluid propagate toward 
the poles on time scales of roughly 500 days. In the Sun, similar patterns are observed in the poleward branch 
of the torsional oscillations, and these may represent a response to poleward propagating magnetic field deep 
below the solar surface. 

Subject headings: convection - MHD - stars:interiors - stars:rotation - stars: magnetic fields - Sun:interior 



1 . STELLAR MAGNETISM AND ROTATION 

Most stars are born rotating quite rapidly. They can ar- 
rive on the main sequence with rot ational velocities as high 
as 200 km s" 1 dBouvier et all 11997b . Stars with convec- 
tion zones at their surfaces, like the Sun, slowly spin down 
as they shed angular momentum t hrough their magn etized 
stellar winds (e.g.. IWeber & Davisl 119671: ISkumanichl 119721: 
iMacGregor & Brenne rl 1 1 99*1F The time needed for signifi- 
cant s pindown appears to be a stron g function of stellar mass 
(e.g.. iBarne s 2003; West et al.ll2.004t) : solar-mass stars slow 
less rapidly than somewhat less massive G and K-type stars, 
but still appear to lose much of their angular momentum by 
the time they are as old as the Sun. Present day observations 
of the solar wind likewise indicate that the current angular 
momentum flux from the Sun is a few times 10 3l) dyn cm (e.g., 
iPizzo et alii 19831) . suggesting a time scale for substantial an- 
gular momentum loss of a few billion years. Thus the Sun 
likely rotated significantly more rapidly in its youth than it 
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does today. 



1.1. Rotation-Activity Relations 



Rotation appears to be inextricably linked to stellar mag- 
netic activity. Observations indicate that in stars with exten- 
sive convective envelopes, chromospheric and coronal activity 
- which partly trace magnetic heating of stellar atmospheres - 
first rise with increasing rotation rate, then eventually level off 
at a constant value for rot ation rates above a mass-dependent 
threshold velocity 



1998; Pizz olato et alJ2003l) . Activit y may eve n decline some- 
what in the most rapid rotators (e.g. JJames et al.l l2000). Sim- 
ilar correspondence is observed between rotation rate and es- 
timates of the unsigned surface magnetic flux (Saarl [l996l 
1200 U iReiners et al.ll2009b . This "rotation-activity" relation- 
ship is tightened when stellar rotation is given in terms of 
the Rossby number Ro ~ P/r c , with P the rotation period 
and t c an estimate of the convective overturning time (e.g., 
iNoves et"ail 119841) . Expressed in this fashion, a common 
rotation-activity correlation appears to span spectral types 



Noves et al.1 11984 iDelfosse et all 
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ranging from late F to late M (e.g.. iMohantv & Basrill2003b 
iPizzolato et aDl2003b iReiners & Basrill2007h . Magnetic fields 
can likewise feed back upon stellar rotation by modifying the 
rate a t which angular momentum is lost through a s tellar wind 
(e.g., IWeber & Davisl Il967t iMatt & Pudritzl 120081) . Analy- 
ses of stellar spindown as a function of age and mass have 
thus provided further constraints on stellar magnetism and its 
connections to rotation. In addition, recent observations of 
solar-type stars may indicate that topology of the global-scale 
fields changes with rotation rate, with the rapid rotators hav- 
ing substan tial global-sca le toroidal magnetic fields at their 
surfaces dPetit et alJl2008l) . The overall picture that emerges 
from these observations is that rapid rotation, as realized in 
the younger Sun and in a host of other stars, can aid in the 
generation of strong magnetic fields and that young stars tend 
to be rapidly rotating and magnetic ally active, whereas older 
ones are slo wer and less active (e.g., Barnes 2003s; IWest et alj 
l2004ll2008h . 

A full theoretical understanding of the rotation-activity re- 
lationship, and likewise of stellar spindown, has remained 
elusive. Some aspects of these phenomena probably de- 
pend upon the details of magnetic flux emergence, chro- 
mospheric and coronal heating, and mass loss mechanisms 
- but the basic existence of a rotation-activity relationship 
is widely thought to reflect some underlyin g rotational de- 
pendence of the dynamo process itself (e.g., iKnobloch et alJ 
[19811 iNoves eTaIlll984t iBaliunas et al.lfl996h . 

1.2. Elements of Global Dynamo Action 

In stars like the Sun, the global-scale dynamo is gener- 
ally thought to be located in the tachocline, an interface of 
shear between the differentially rotating convection zone and 
the radiative interior which is in solid body rotation (e.g., 
Parkerlll993l ICharbonneau & MacGregodl 19971: lOssendrijverl 
2003b . Helioseismology revealed the internal rotation pro- 
file o f the Sun and the prese nce of this important shear layer 
(e.g., [Thompson et al. 2003). The stably stratified tachocline 
may also provide a region for storing and amplifying coherent 
tubes of magnetic field which may eventually rise to the sur- 
face of the Sun as sunspots. It has generally been believed that 
magnetic buoyancy instabilities may prevent fields from being 
strong ly amplified within the bulk of the convection zone it- 
self (lParkeri(l975l) . In the now prevalent "interface dynamo" 
model, solar magnetic fields are partly generated in the con- 
vection zone by helical convection, then transported down- 
ward into the tachocline where they are organized and ampli- 
fied by the shear. Ultimately the fields may become unstable 
and rise to the surface. 

Although the rotational dependence of this process is not 
well understood, some guidance may come from mean-field 
dynamo theory. In such theories, the solar dynamo is re- 
ferred to as an "a— f2" dynamo, with the a-effect char acteriz- 
ing the twisting of fiel ds by helical convection (e.g., Moffatt 
1 19781: ISteenbeck et al.|[l966T) . and the O-effect representing 
the shearing of poloidal fields by differential rotation to form 
toroidal fields. Both of these effects are, in mean-field theory, 
sensitive to rotation - the a-effect because it is proportional 
to the kinetic helicity of the convective flows, which sense the 
overall rotation rate, and the £l-effect because more rapidly 
rotating stars are generally expected to have stronger differ- 
ential rotation. But the detailed nature of these effects in the 
solar dynamo and the appropriate scaling with rotation has 
been very difficult to elucidate. 

Simulations of the global-scale solar dynamo have gener- 



ally affirmed the view that the tachocline may play a cen- 
tral role in building the globally-ordered magnetism in the 
Sun. Early three-dimensional (3D) simulations of solar con- 
vection without a tachocline at the base of the convection zone 
achieved dynamo action and produced magnetic fields which 
were strongly dominat ed by fluctuating components with lit- 
tle global-scale order dBrun et alj 12004). When a tachocline 
of penetration and shear was included, remarkable global- 
scale structures were realized in the tachocline region, while 
the convection zone r emained dominated by fluctuating fields 
dBrowning et alj|2006l) . These simulations are making good 
progress toward clarifying the elements at work in the opera- 
tion of the solar global-scale dynamo, but for other stars many 
questions remain. In particular, observ ations of large-scal e 
magnetism in fully convective M-stars (Do nati et al.1 2006), 
along with the persistence of a rotation-activity correlation in 
such low-mass stars, hint that perhaps tachoclines may not be 
essential for the generation of global-scale magnetic fields. 
This view is partly borne out by simulations of M -dwarfs 
under strong rotational constraints (iBrowningl l2008h . where 
strong longitudinal mean fields were realized despite the lack 
of either substantial differential rotation or a stable interior 
and thus no classical tachocline. Major puzzles remain in 
the quest to understand stellar magnetism and its scaling with 
stellar rotation. 

1.3. Convection and Dynamos in Rapidly Rotating Systems 

We began our study of rapidly rotating suns by carrying 
out a suite of 3D hydrodynamic simulations in full spherical 
shells that explored the coupling o f rotation and co nvection 
in these younger solar-type stars dBrown et alJl2008h . Those 
simulations studied the influence of rotation on the patterns of 
convection and the nature of global-scale flows in such stars. 
The shearing flows of differential rotation generally grow in 
amplitude with more rapid rotation, possessing rapid equators 
and slower poles, while the meridional circulations weaken 
and break up into multiple cells in radius and latitude. More 
rapid rotation can also substantially modify the patterns of 
convection in a surprising fashion. With more rapid rotation, 
localized states begin to appear in which the convection at 
low latitudes is modulated in its strength with longitude. At 
the highest rotation rates, the convection can become confined 
to active nests which propagate at distinct rates and persist for 
long epochs. 

Motivated by these discoveries, we turn here to explorations 
of the possible dynamo action achieved in solar-type stars ro- 
tating at three and five times the current solar rate. These 3D 
magnetohydrodynamic (MHD) simulations span the convec- 
tion zone alone, as the nature of tachoclines in more rapidly 
rotating suns is at present unclear. We shall show that a va- 
riety of dynamos can be excited, including steady and oscil- 
lating states, and that dynamo action is substantially easier to 
achieve at these faster rotation rates than in the solar simu- 
lations. Magnetism leads to strong feedbacks on the flows, 
particularly modifying the differential rotation and its scal- 
ing with overall rotation rate Slo- The magnetic fields which 
form in these dynamos have prominent global-scale organiza- 
tion within the conv ection zone, in contrast to previous sola r 
dynamo simulations (iBrun et alJ2004t [Browning et al.ll2006l) . 
Quite strikingly, we find that coherent global magnetic struc- 
tures arise naturally in the midst of the turbulent convection 
zones. These wreath-like structures are regions of strong lon- 
gitudinal field Z?0 organized loosely into tubes, with fields 
wandering in and out of the surrounding convection. These 
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wreathes of magnetism differ substantially from the idealized 
flux tubes supposed in many dynamo theories, though they 
may be related to coherent structures achi eved in local sim- 
ulations of dynamo action in s hear flows dCline et al J 120031: 
IVasil & Bru mmell 2008, 2009). Here we explore the nature 
of magnetic wreathes realized in our global simulations, and 
discuss their temporal behavior. 

We outline in $2] the 3D MHD anelastic spherical shell 
model and the parameter space explored by these simula- 
tions. We then examine in §^3]and|4]the structure of magnetic 
fields found in rapidly rotating dynamo at three times the so- 
lar rate, which builds remarkable global-scale ordered fields 
in the midst of its convection zone. In §§0|6]we explore time 
varying behavior and organized global-scale polarity rever- 
sals in a dynamo rotating at five times the current solar rate. 
In §|7]we return to the three solar case and examine how such 
global-scale fields are created and maintained. We reflect on 
our findings in ^8] 

2. GLOBAL MODELLING APPROACH 

To study the coupling between rotation, magnetism and 
the large-scale flows achieved in stellar convection zones, we 
must employ a global model which simultaneously captures 
the spherical shell geometry and admits the possibility of 
zonal jets and large eddy vortices, and of convective plumes 
that may span the depth of the convection zone. The solar 
convection zone is intensely turbulent and microscopic values 
of viscosity and magnetic and thermal diffusivities in the Sun 
are estimated to be very small. Numerical simulations cannot 
hope to resolve all scales of motion present in real stellar con- 
vection and must instead strike a compromise between resolv- 
ing dynamics on small scales and capturing the connectivity 
and geometry of the global scales. Here we focus on the latter 
by studying a full spherical shell of convection. 

2.1. Anelastic MHD Formulation 



Our tool for exploring MHD stellar convection is the 
anelastic s pherical harmonic ( ASH) code, which is described 
in detail in lClune et ap (119991). The imp li 
Bran et all (12004). 



netism is discussed in 



lementation of mag- 
ASH solves the 3D 



MHD anelastic equations of motion in a rotating spherical 
shell using the pseudo-spectral method and runs efficiently on 
massively parallel architectures. We use the anelastic approx- 
imation to capture the effects of density stratification with- 
out having to resolve sound waves which have short periods 
(about 5 minutes) relative to the dynamical time scales of the 
global scale convection (weeks to months) or possible cycles 
of stellar activity (years to decades). This criteria effectively 
filters out the fast magneto-acoustic modes while retaining the 
slow modes and Alfven waves. Under the anelastic approxi- 
mation the thermodynamic fluctuating variables are linearized 
about their spherically symmetric and evolving mean state, 
with radially varying density p, pressure P, temperature f and 
specific entropy S. The fluctuations about this mean state are 
denoted as p, P, T and S. In the reference frame of the star, 
rotating at average rotation rate f^o, the resulting MHD equa- 
tions are: 

V-(pv) = 0, (1) 

(2) 



dv 
dt 



VB = 0, 
+ (v- V)v + 2fi X v 



-V(P+P) 



1 



(3) 



0B 

— = V x (v x B)-Vx (77V x B), 
at 



pf 



OS 

— + v-V(S+S) 



V- [K r pc p V(T + T) + K Q pTVS+KpTVS] 



(4) 



(5) 



f + 2pv 



njeij--{V-vY 



where v = (v r ,V6,v#) is the local velocity in the stellar ref- 
erence frame, B = (B r ,Be,B,p) is the magnetic field, j is the 
vector current density, g is the gravitational acceleration, c p is 
the specific heat at constant pressure, n r is the radiative diffu- 
sivity and T> is the viscous stress tensor, given by 



T>ij = —2pv 



eij--{V-v)8ij 



(6) 



where ey is the strain rate tensor. Here v, k and rj are the dif- 
fusivities for vorticity, entropy and magnetic field. We assume 
an ideal gas law 

P = KpT, (7) 

where 1Z is the gas constant, and close this set of equations 
using the linearized relations for the thermodynamic fluctua- 
tions of 

p _P T _ P S 
~p~P~T ~jP~~ p ' 



(8) 



The mean state thermodynamic variables that vary with radius 
are evolved with the fluctuations, thus allowing the convection 
to modify the entropy gradients which drive it. 

The mass flux and the magnetic field are represented with a 
toroidal-poloidal decomposition as 

pv = Vx Vx(Wr) + Vx(Zr), (9) 
B= Vx Vx(/3r) + Vx(Cr), (10) 

with streamfunctions W and Z and magnetic potentials j3 
and This approach ensures that both quantities remain 
divergence-free to machine precision throughout the simula- 
tion. The velocity, magnetic and thermodynamic variables 
are all expanded in spherical harmonics for their horizon- 
tal structure and in Chebyshev polynomials for their radial 
structure. The solution is time evolved with a second-order 
Adams-Bashforth/Crank-Nicolson technique. 

ASH is a large-eddy simulation (LES) code, with subgrid- 
scale (SGS) treatments for scales of motion which fall be- 
low the spatial resolution in our simulations. We treat these 
scales with effective eddy diffusivities, v, k and 77, which 
represent the transport of momentum, entropy and magnetic 
field by unresolved motions in the simulations. These sim- 
ulations are based on the hydrodynamic studies reported in 
iBrown et al.l (|2008), and as there v, k and rj are taken for 
simplicity as functions of radius alone and proportional to 
p~ x l\ This adopted SGS variation, as in iBrun et all (120041) 
and lBrowning et all d2006l) . yields lower diffusivities near the 
bottom of the layer and thus higher Reynolds numbers. Act- 
ing on the mean entropy gradient is the eddy thermal diffusion 
/to which is treated separately and occupies a narrow region 
in the upper convection zone. Its purpose is to transport en- 
tropy through the outer surface where radial convective mo- 
tions vanish. 



+(P+P)9+ -;- ( V X B) X B- V ■ V, 
Air 
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TABLE 1 

Parameters for Primary Simulations 



Case N r ,N e ,N^ Ra Ta Re Re' Rm Rm' Ro Roc v r) n /fl e 

D3 96x256x512 3.22xl0 5 1.22xl0 7 173 105 86 52 0.378 0.311 1.32 2.64 3 

D5 96x256x512 1.05 xlO 6 6.70xl0 7 273 133 136 66 0.273 0.241 0.940 1.88 5 

H3 96x256x512 4.10xl0 5 1.22xl0 7 3 35 1 05 — — 0.427 0.3 5 3 1.32 — 3 

H5 96x256x512 1.27xl0 6 6.70xl0 7 576 141 — — 0.303 0.268 0.940 — 5 



NOTE. — Dynamo simulations at three and five times the solar rotation rate are cases D3 and D5, and their hydrodynamic 
(non-magnetic) companions are H3 and H5. All simulations have inner radius rt, ot = 5.0 x 10 I0 cm and outer radius of r top = 
6.72 x 10 lo cm, with L = (rt p-n, ot ) = 1.72 X 10 I0 cm the thickness of the spherical shell. Evaluated at mid-depth are the 
Rayleigh number Ra= (-dp/dS)(dS/dr)gL*/pun, the Taylor number Ta = 4f! 2 ,L 4 /i/ 2 , the rms Reynolds number Re = v !ras L/v 
and fluctuating Reynolds number Re' = v' Yms L/v, the magnetic Reynolds number Rm = VnnsL/r) and fluctuating magnetic 
Reynolds number Rm' = v' rms L/r], the Rossby number Ro = lo/2Qq , and the convective Rossby number Roc = (Ra/TaPr) 1 / 2 . 
Here the fluctuating velocity v' has the axisymmetric component removed: v' = v— (v), with angle brackets denoting an average 
in longitude. For all simulations, the Prandtl number Pr = v/k is 0.25 and the magnetic Prandtl number Pm = v/r) is 0.5. The 
viscous and magnetic diffusivity, v and r\, are quoted at mid-depth (in units of 10 12 cm 2 s~'). The rotation rate f!o of each 
reference frame is in multiples of the solar rate Qq = 2.6 X 10~ 6 rads~' or 414 nHz. The viscous time scale at mid-depth 
t v = L 2 jv is 3640 days for case D5 and the resistive time scale is 1820 days. Rotation periods at three and five times the solar 
rate are in turn 9.3 days and 5.6 days. 



The boundary conditions imposed at the top and bottom of 
the convective unstable shell are: 

1 . Impenetrable top and bottom: v r = , 

2. Stress-free top and bottom: 

(d/dr)(v e /r) = (d/dr)(v^/r) = 0, 

3. Constant entropy gradient at top and bottom: 

d(S+ S)/dr= const, (11) 

4. Match to external potential field at top: 

B = V$ and V 2 $ = 0|, ,. , 

1 flop : 

5. Perfect conductor at bottom: 

B r = (d/dr)(B e /r) = (fl/flr)( V) = . 

2.2. Posing the Dynamo Problem 

Our simulations are a simplified picture of the vastly tur- 
bulent stellar convection zones present in G-type stars. We 
take solar values for the input entropy flux, mass and radius, 
and explore simulations of two stars rotating at three and five 
times the current solar rotation rate. We focus here on the 
bulk of the convection zone, with our computational domain 
extending from 0.72/? Q to 0.97/?q, thus spanning 172 Mm in 
radius. The total density contrast across the shell is about 25. 
The reference or mean state of our thermod ynamic variables 
is derived from a ID solar structure model dBrun et alj|2002l) 
and is continuously updated with the spherically symmetric 
components of the thermodynamic fluctuations as the simula- 
tions proceed. The refe rence state in all of t hese simulations is 
similar to that shown in Brown et al. (2008). We avoid regions 
near the stellar surface where hydrogen ionization and radia- 
tive losses drive intense convection (like granulation) on very 
small scales that we cannot resolve, and thus position the up- 
per boundary slightly below this region. Our lower boundary 
is positioned near the base of the convection zone, thus omit- 
ting the stably stratified radiative interior and the shear layer 
at the base of the convection zone known as the tachocline. 



As in our previous work (iBrown et alJ l2008). we conducted 
our simulations at three and five times the solar rate f2© on a 

—2/3 

path where the SGS diffusivities v, k and 77 decrease as fl Q , 
in order to maintain vigorous convection as rotation attempts 
to constrain and quench the motions. The fundamental char- 
acteristics of our simulations and parameters definitions are 
summarized in Table Q] 

The dynamo simulations were initiated from mature hy- 
drodynamic progenitor cases which had been evolved for at 
least 5000 days at each rotation rate and were well equili- 
br ated. The progenitor s are very similar to the cases reported 
in IBrown et alj d2008l) . but we chose a functional form for 
the SGS entropy diffusion kq that is more confined to the 
upper 10% of the convection zone; the unresolved flux here 
does not vary with rotation rate. The effects of this are sub- 
tle, resulting primarily in slightly stronger latitudinal gradi- 
ents of differential rotation and temperature in the uppermost 
regions of the shell. The patterns of convection are very sim- 
ilar, though slightly more complex near the top of the shell, 
and the Reynolds number remains high throughout the shell. 
These cases were well evolved and possess intricate convec- 
tive patterns and solar-like differential rotation profiles, with 
fast zonal flow at the equator and slower flows at the poles. 

To initiate our dynamo cases, a small seed dipole magnetic 
field was introduced and evolved via the induction equation. 
The energy in the magnetic fields is initially many orders of 
magnitude smaller than the energy contained in the convective 
motions, but these fields are amplified by shear and grow to 
become comparable in energy to the convective motions. 

These dynamo simulations are computationally intensive, 
requiring both high resolutions to correctly represent the ve- 
locity fields and long time evolution to capture the equili- 
brated dynamo behavior, which may include cyclic variations 
on time scales of several years. The strong magnetic fields can 
produce rapidly moving Alfven waves which seriously restrict 
the Courant-Friedrichs-Lewy (CFL) timestep limits in the up- 
per portions of the convection zone. Case D3, rotating three 
times faster than the current Sun, has been evolved for over 
7000 days (or over 2 million timesteps), and case D5, rotating 
five times faster than the Sun, has seen more than 17000 days 
of evolution (representing more than 10 million timesteps). 
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FIG. 1. — Corrective structures and mean flows in cases D3 and H3. (a) Radial velocity v, in dynamo case D3, shown in global Mollweide projection at 
O.95R0, with upflows light and downflows dark. Poles are at top and bottom and the equator is the thick dashed line. The stellar surface at Rq is indicated by 
the thin surrounding line, (b) Profiles of mean angular velocity Q(r,8), accompanied in (c) by radial cuts of Q at selected latitudes. A strong differential rotation 
is established by the convection, (d) Profiles of meridional circulation, with sense of circulation indicated by color (red counter-clockwise, blue clockwise) and 
streamlines of mass flux overlaid, (e-h) Companion presentation of fields for hydrodynamic progenitor case H3. The patterns of radial velocity are very similar 
in both cases. The differential rotation is much stronger in the hydrodynamic case and the meridional circulations there are somewhat weaker, though their 
structure remains similar. 



We plan to report on a variety of other dynamo cases, some 
at higher turbulence levels and rotation rates, in a subsequent 
paper. 

These two cases were conducted at magnetic Prandtl num- 
ber Pm =v jr\ = 0.5, a value significantly lower tha n employed 
in our previous solar simulations. I n particular. iBrun et alj 
( 2004) explored Pm = 2 , 2 .5 and 4, and lBrowning et alj d2006l) 
studied Pm = 8. The high magnetic Prandtl numbers were re- 
quired in the solar simulations to reach sufficiently high mag- 
netic Reynolds num bers to drive susta ined dynamo action. 
In the simulations of iBrun et al.1 (|2004) only the simulations 
with Pm > 2.5 and Rm' > 300 achieved sustained dynamo ac- 
tion, where Rm' is the fluctuating magnetic Reynolds number. 
We are here able to use a lower magnetic Prandtl number for 
three reasons. Firstly, more rapid rotation tends to stabilize 
convection and lower values of v and 77 are required to drive 
the convection. Once convective motions begin however they 
become quite vigorous and the fluctuating velocities saturate 
at values comparable to our solar cases. Thus the Reynolds 
numbers achieved are fairly large and we can achieve mod- 
estly high magnetic Reynolds numbers even at low Pm. Sec- 
ondly, the differential rotation becomes substantially stronger 
with both more rapid rotation f2n and with lower diffusivities 
v and 77. This global-scale flow is an important ingredient 
and reservoir of energy for these dynamos, and the increase 
in its amplitude means that low Pm dynamos can still achieve 
large magnetic Reynolds numbers based on this zonal flow. 
Lastly, the critical magnetic Reynolds number for dynamo 
action likely decr eases with increasing kinetic helicity (e.g., 
Leora t etld1ll981l). Helicity g enerally increases with rotation 
rate (e.g.. iKapvla et al.l 12009). so the rapidly rotating flows 
considered here achieve dynamo action at somewhat lower 



Rm than the models of IBrun et al.l (120041) . which rotated at 
the solar rate. 

3. DYNAMOS WITH PERSISTENT MAGNETIC WREATHES 

We begin by turning to case D3 which yielded fairly persis- 
tent wreathes of magnetism in its two hemispheres, though 
these did wax and wane somewhat in strength once estab- 
lished. Examining the properties of this dynamo solution help 
to provide a perspective for the greater variations realized in 
our more rapidly rotating case D5. 

3.1. Patterns of Convection 

The complex and evolving convective structures in our dy- 
namo cases are substantially similar to the patterns of con- 
vection found in our hydrodynamic simulations. Our dynamo 
solution rotating at three times the solar rate, case D3, is pre- 
sented in Figure Q] along with its hydrodynamic progenitor, 
case H3. The radial velocities shown near the top of the sim- 
ulated domain (Figs. Qj!,e) have broad upflows and narrow 
downflows as a consequence of the compressible motions. 
Near the equator the convection is aligned largely in the north- 
south direction, and these broad fronts sweep through the do- 
main in a prograde fashion. The strongest downflows pene- 
trate to the bottom of the convection zone; the weaker flows 
are partially truncated by the strong zonal flows of differential 
rotation. In the polar regions the convection is more isotropic 
and cyclonic. There the networks of downflow lanes surround 
upflows and both propagate in a retrograde fashion. 

The convection establishes a prominent differential rota- 
tion profile by redistributing angular momentum and entropy, 
building gradients in latitude of angular velocity and temper- 
ature. Figures []J>,/ show the mean angular velocity Sl(r, 9) 
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FIG. 2. — Magnetic wreathes and convective flows sampled at the same instant in case D3. (a) Longitudinal magnetic field near the top of the shell (O.95R0) 
and (b) at mid-depth (O.85R0). Strong flux structures with opposite polarity lie above and below the equator and span the convection zone. (c,d) Weaker radial 
magnetic field B r permeates and encircles each wreath. (e,f) Strong convective upflows and downflows shown by V r pass through and around the wreathes. The 
regions of strong magnetism tend to disrupt the convective flows while the strongest downflows serve to pump the wreathes to greater depths. 



for cases D3 and H3, revealing a solar-like structure with a 
prograde (fast) equator and retrograde (slow) pole. This prop- 
erty is also realized for cases D5 and H5 with faster rotation. 
Figures [TJ;,^ present in turn radial cuts of fl at selected lat- 
itudes, which are useful as we consider the angular velocity 
patterns realized here with faster rotation. These fl(r, 9) pro- 
files are averaged in azimuth (longitude) and time over a pe- 
riod of roughly 200 days. Contours of constant angular veloc- 
ity are aligned nearly on cylinders, influenced by the Taylor- 
Proudman theorem. 

In the Sun, helioseismology has revealed that the contours 
of angular velocity are aligned almost on radial lines rather 
than on cylinders. The tilt of il contours in the Sun may be due 
in part to the thermal structure of th e solar tachocli ne, as first 
found in the mean-field models of Rempel (2005) an d then 
in 3D simulations of global-scale convection by Miesc h~et alJ 
(2006). In those computations, it was realized that introduc- 
ing a weak latitudinal gradient of entropy at the base of the 
convection zone, consistent with a thermal wind balance in 
a tachocline of shear, can serve to tilt the £1 contours toward 
a more radial alignment without significantly changing either 
the overall £1 contrast with latitude or the convective patterns. 
We expect similar behavior here, but at present, observations 
of rapidly rotating stars only measure differential rotation at 



the surface and do not offer constraints on either the exis- 
tence of tachoclines in young suns or the nature of their inter- 
nal differential rotation profiles. As such, we have neglected 
the possible tachoclines of penetration and shear entirely in 
these models and instead adopt the simplification of imposing 
a constant radial entropy gradient at the bottom of the convec- 
tion zone. 

The differential rotation achieved is stronger in our hydro- 
dynamic case H3 than in our dynamo case D3. This can be 
quantified by measurem ents of the lat i tudina l angular velocity 
shear Afii at . Here, as in lBrown et al.l J2008), we define Afii at 
as the shear near the surface between the equator and a high 
latitude, say ±60° 

Afii a t = (12) 

and the radial shear ASl r as the angular velocity shear between 
the surface and bottom of the convection zone near the equator 

Afi r = f2o.97fi - ^O.72R • (13) 

We further define the relative shear as AS7i at /S7 eq . In both 
definitions, we average the measurements of Aft in the north- 
ern and southern hemispheres, as the rotation profile is often 
slightly asymmetric about the equator. Case H3 achieves an 
absolute contrast Afii at of 2.22 /1 rads" 1 (352 nHz) and a rela- 
tive contrast of 0.247. The strong global-scale magnetic fields 
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TABLE 2 
Near-surface AO 



FIG. 3. — Field line tracings of magnetic wreathes in case D3. (a) Snapshot of two wreathes in full volume at same instant as in Fig. [2] Lines trace the magnetic 
fields, color denoting the amplitude and polarity of the longitudinal field (red, positive; blue, negative). Magnetic field threads in and out of the wreathes, 
connecting the two opposite polarity structures across the equator (i.e., region A) and to the polar regions where the magnetic field is wound up by the cyclonic 
convection, (b) Same snapshot showing south polar region, (c) Zoom in on region A showing the complex interconnections across the equator between the two 
wreathes and to high latitudes. Convective flows create the distinctive waviness visible in all three images. 

3.2. Kinetic and Magnetic Energies 

Convection in these rapidly rotating dynamos is responsible 
for building the differential rotation and the magnetic fields. 
In a volume averaged sense, the energy contained in the mag- 
netic fields in case D3 is about 10% of the kinetic energy. 
About 35% of this kinetic energy is contained in the fluc- 
tuating convection (CKE) and about 65% in the differential 
rotation (DRKE), whereas the weaker meridional circulations 
contain only a small portion (MCKE). The magnetic energy is 
split between the contributions from fluctuating fields (FME), 
involving roughly 53% of the total magnetic energy, and the 
energy of the mean toroidal fields (TME) that are 43% of the 
total. The energy contained in the mean poloidal fields (PME) 
is only 4% of the total magnetic energy. These energies are 
defined as 



Case 


Afii at 


AO r 


Afiiat/neq 


Epoch 


D3 


1.17 


0.71 


0.136 


6460-6920 


D5<w« 


1.14 


0.71 


0.083 


3500-5700 


D5""" 


0.91 


0.39 


0.067 


3702 


D5 max 


1.43 


0.98 


0.102 


4060 


H3 


2.22 


0.94 


0.246 




H5 


2.77 


1.31 


0.192 





NOTE. — Angular velocity shear in units of 
/irads -1 , with Af2i al and Af2i at /f2 C q measured near 
the surface (0.97Sq) and AC, measured across the 
full shell at the equator. For the dynamo cases, these 
measurements are taken over the indicated range of 
days. In oscillating case D5, these measurements 
are averaged over a long epoch (avg), and are also 
taken at two instants in time when the differential ro- 
tation is particularly strong (max) and when magnetic 
fields have suppressed this flow (miri). The hydrody- 
namic cases are each averaged for roughly 300 days. 
Case D3 also shows slow variations in Af2i at over pe- 
riods of about 2000 days. 



realized in the dynamo case D3 serve to diminish the differen- 
tial rotation. As such, this case achieves an absolute contrast 
Af2i at of only 1.17 n rads" 1 (186 nHz) and a relative contrast 
of 0.136. This results from both a slowing of the equatorial 
rotation rate and an increase in the rotation rate in the polar re- 
gions. These results are quoted in Table|2] along with related 
measurements for our five solar dynamo case D5 and hydro- 
dynamic cases H3 and H5. It is interesting to note that in 
both dynamo cases the amount of latitudinal and radial shear 
is almost the same, whereas in the hydrodynamic simulations 
the more rapidly rotating case H5 has larger angular velocity 
contrasts than the slower case H3. 

The meridional circulations realized in the dynamo case D3 
are very similar to those found in its hydrodynamic progeni- 
tor (case H3). As illustrated in Figures ld,h, the circulations 
are multi-celled in radius and latitude. The cells are strongly 
aligned with the rotation axis, though some flows along the in- 
ner and outer boundaries cross the tangent cylinder and serve 
to couple the polar regions to the equatorial convection. Flows 
of meridional circulation are slightly stronger in the dynamo 
cases than in the purely hydrodynamic cases, and these flows 
weaken with more rapid rotation. 



CKE = 



1 _ 

2 P 



(v r -(Vr)) +(v e - 



M 2 + 



DRKE=-p(v ) 2 , 
MCKE=ip((v r ) 2 +(v e ) 2 



FME= — 



1 

8^ 



(B r -(B r )) 2 +(Bg-(B e )Y + 



{B*-(B*)Y 



TME: 



PME = 



■sf(*> 



■+06 



(14) 
(15) 
(16) 

(17) 
(18) 
(19) 



where angle brackets denote an average in longitude. 

These results are in contrast to our previous simulations of 
the solar dynamo, where the mean fields contained only about 
2% of the m agnetic energy an d the fluctuating fields contained 
nearly 98% dBrun et al. 2004). In simulations of the solar dy- 
namo that inc luded a stable tachocli ne at the base of the con- 
vection zone dBrowning et al.i r2006). the energy of the mean 
fields in the tachocline can exceed the energy of the fluctuat- 
ing fields there by about a factor of three, though the fluctuat- 
ing fields still dominate the magnetic energy budget within the 
convection zone itself. Simulations of dynamo act ivity in the 
convecting cores of A-type stars dBrun et al.ll2005l) achieved 
similar results. There in the stable radiative zone the energies 
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TABLE 3 
Energies 



Fig. 4. — Persistent wreathes of magnetism in case D3. (a) Time-latitude plots of azimuthally-averaged longitudinal field (Bj,) at mid-convection zone 
(O.85/?0) in a view spanning latitudes from ±70°, with scaling values indicated. The two wreathes of opposite polarity persist for more than 4000 days. 
(b) Mean colatitudinal magnetic field (Bg) at mid-convection zone over same interval, (c) Snapshots of B^ in Mollweide projection at mid-convection zone, 
shown for three times indicated in a.b. The wreathes maintain constant polarity over long time intervals, but still show variation as they interact with the 
convection. Time t2 corresponds to the snapshot in Fig.f2j). 

very different rotation rates. This is in striking contrast to the 
behavior of the hydrodynamic cases, where DRKE is much 
larger in the more rapidly rotating case H5 than in case H3. 
Meridional circulations are comparably weak in all cases. 

The amount of energy contained in the magnetic fields is 
different in these two dynamo cases, with energies generally 
stronger in case D5 than in case D3. In an average sense, all 
three magnetic energies are about 1.5 times greater in case D5 
than in case D3. Case D5 shows substantial time variation, 
and at periods either FME or TME can become quite simi- 
lar to those values realized in case D3. Meanwhile, PME is 
always stronger in D5 than in D3. 



Case 


CKE 


DRKE 


MCKE 


FME 


TME 


PME 


D3 


2.31 


4.35 


0.010 


0.36 


0.29 


0.029 


D5 m '« 


1.85 


4.46 


0.006 


0.55 


0.43 


0.048 


D5""" 


1.70 


2.85 


0.005 


0.50 


0.25 


0.062 


D5 max 


1.85 


7.52 


0.007 


0.39 


0.65 


0.042 


H3 


2.56 


22.2 


0.012 








H5 


2.27 


34.3 


0.008 









NOTE. — Volume-averaged energy densities relative to the 
rotating coordinate system. Kinetic energies are shown for 
convection (CKE), differential rotation (DRKE) and merid- 
ional circulations (MCKE). Magnetic energies are shown for 
fluctuating magnetic fields (FME), mean toroidal fields (TME) 
and mean poloidal fields (PME). All energy densities are re- 
ported in units of 10°ergcm~ 3 and are averaged over 1000 day 
periods except for time-varying case D5, where intervals are 
as in Tablef2] 



of the mean fields were able to exceed the energy contained 
in the fluctuating fields, but in the convecting core the fluc- 
tuating fields contained roughly 95% of the magnetic energy. 
Simulations of dynamo action in fully-convective M-stars do 
howev er show high lev els of magnetic energy in the mean 
fields (Browning 2008). In those simulations the fluctuating 
fields still contain much of the magnetic energy, but the mean 
toroidal fields possess about 18% of the total throughout most 
of the stellar interior. In our rapidly rotating suns, the mean 
fields comprise a significant portion of the magnetic energy 
in the convection zone and are as important as the fluctuating 
fields. 

Convection is similarly strong in all four rapidly rotating 
cases, and CKE is similar in magnitude. The differential ro- 
tation in the dynamo cases is much weaker than in the hy- 
drodynamic progenitors and DRKE has decreased by about a 
factor of five. The magnetic fields and differential rotation in 
case D5 change in time, but the average energy contained in 
DRKE is nearly the same in both dynamo cases despite their 



4. WREATHES OF MAGNETISM 

These dynamos produce striking magnetic structures in the 
midst of their turbulent convection zones. The magnetic field 
is organized into large banded, wreath-like structures posi- 
tioned near the equator and spanning the depth of the con- 
vection zone. These wreathes are shown at two depths in 
the convection zone in Figure |2] The dominant component 
of the magnetic wreathes is the strong longitudinal field B<p, 
with each wreath possessing its own polarity. The average 
strength of the longitudinal field at mid-convection zone is 
±7 kG and peak field strengths there reach roughly ±26 kG. 
Threaded throughout the wreathes are weaker radial and lat- 
itudinal magnetic fields, which connect the two structures 
across the equator and also to the high-latitude regions. 

These wreathes of magnetism survive despite being em- 
bedded in vigorous convective upflows and downflows. The 
convective flows leave their imprint on the magnetic struc- 
tures, with individual downflow lanes entraining the magnetic 
field, advecting it away, and stretching it into B r while leav- 
ing regions of locally reduced B^. The slower upflows carry 
stronger B^ up from the depths. Where the magnetism is 
particularly strong the convective flows are disrupted. Mean- 
while, where the convective flows are strongest, the longitu- 
dinal magnetic field is weakened and appears to vanish. In re- 
ality, the magnetic wreathes here are diving deeper below the 
mid-convection zone, apparently pumped down by the pum- 
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meling action of the strong downfiows. 

The deep structure of these wreathes is revealed by field 
line tracings throughout the volume, shown in Figure [3] for 
the same instant in time. The wreathes are topologically leaky 
structures, with magnetic field lines threading in and out of the 
surrounding convection. The wreathes are connected to the 
high-latitude (polar) convection, and on the poleward edges 
they show substantial winding from the highly vortical con- 
vection found there. This occurs in both the northern and 
southern hemispheres, as shown in two views at the same in- 
stant (north, Fig. [3j3t and south, Fig. |3j>). It is here that the 
global-scale poloidal field is being regenerated by the cou- 
pling of fluctuating velocities and fluctuating fields. Magnetic 
fields cross the equator, tying the two wreathes together at 
many locations (Fig. [3]:). The strongest convective downfiows 
leave their imprint on the wreathes as regions where the field 
lines are dragged down deeper into the convection zone, yield- 
ing a wavy appearance to the wreathes as a whole. 

4. 1 . Wreathes Persist for Long Epochs 

The wreathes of magnetism built in case D3 persist for long 
periods of time, with little change in strength and no rever- 
sals in global-scale polarity for as long as we have pursued 
these calculations. The long-term stability of the wreathes re- 
alized by the dynamo of case D3 is shown in Figure |4] Here 
the azimuthally-averaged longitudinal field (B$) and colatitu- 
dinal field (Bg) are shown at mid-convection zone at a point 
after the dynamo has equilibrated and for a period of roughly 
5000 days (i.e., several ohmic diffusion times). During this in- 
terval there is little change in either the amplitude or structure 
of the mean fields. This is despite the short overturn times of 
the convection (10-30 days) or the rotation period of the star 
(~ 9 days). The ohmic diffusion time at mid-convection zone 
is approximately 1300 days. 

Though the mean (global-scale) fields are roughly steady in 
nature (Figs.lU^b), the magnetic field interacts strongly with 
the convection on smaller scales. Several samples of longitu- 
dinal field Bff, are shown in full Mollweide projection at mid- 
convection zone (Fig.HJ;). The magnetic fields are clearly re- 
acting on short time scales to the convection but the wreathes 
maintain their coherence. 

There are also some small but repeated variations in the 
global-scale magnetic fields. Visible in Figure |4}? are events 
where propagating structures of (Bg) reach toward higher lat- 
itudes over periods of about 1000 days (i.e., from day 3700 
to day 4500 and from day 5600 to day 6400). These are ac- 
companied by slight variations in the volume-averaged mag- 
netic energy densities and the comparable kinetic energy of 
the differential rotation. These variations are also visible in 
the differential rotation itself, as shown in Figure The dif- 
ferential rotation is fairly stable, though some time variation 
is visible at high latitudes. This is better revealed (Fig.[5J>) by 
subtracting the time-averaged profile of Q, at each latitude, re- 
vealing the temporal variations about this mean. In the polar 
regions above ±40° latitude, speedup features move poleward 
over 500 day periods. These features track similar structures 
visible in the mean magnetic fields (Fig.@J?). 

These evolving structures of magnetism and faster and 
slower differential rotation appear to be the first indications 
of behavior where the mean fields themselves begin to wax 
and wane substantially in strength. As the magnetic Reynolds 
number is increased, this time varying behavior becomes 
more prominent and can even result in organized changes 
in the global-scale polarity. Such behavior is evident in our 
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FIG. 5. — Differential rotation in case D3. (a) Angular velocity f! at mid- 
convection zone (0.85i?o), with ranges in both nHz and /xrads -1 . The equator 
is fast while the poles rotate more slowly, (b) Temporal variations are empha- 
sized by subtracting the time-averaged profile of Q(r,8), revealing speedup 
structures at high latitudes and pulses of f ast a nd slow motion near the equa- 
tor, (c) Angular velocity shear AQi at (eg. 112) near the surface (upper curve, 
green) and at mid-convection zone (lower, red). 

case D5. 

5. A CYCLIC DYNAMO 

Our more rapidly rotating simulation, case D5 at five times 
the current solar rate, also builds strong wreathes of mag- 
netism that span the convection zone. In this simulation the 
dynamo is not steady in time and instead goes through cycles 
of activity. During a cycle, the global-scale magnetic fields 
wax and wane in strength and at the lowest field strengths 
they can flip their polarity. We shall begin by looking at the 
general properties of the convective flows and their associated 
differential rotation and meridional circulation, and then turn 
to examining the nature of the magnetic fields and their time- 
varying behavior. 

5.1. Patterns of Convection in Case D5 

Figure |6j? shows a snapshot of the patterns of convection 
realized in case D5 near the top of the domain. Much as in 
the radial velocity patterns shown for case D3 (Figs. [Bi|2^), 
here with faster rotation we continue to have prominent north- 
south aligned cells in the lower latitudes and more isotropic 
patterns near the poles. There is some modulation with longi- 
tude in the equatorial roll amplitudes. The downfiows are fast 
and narrow, while the upflows are broader and slower. The 
convection establishes a prominent differential rotation, with 
a fast equator and slow poles (Fig. [(J). 

As in case D3, here too the convective downflow structures 
propagate more rapidly than the differential rotation in which 
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FIG . 6 . — Patterns of convection at five times the solar rotation rate, (a) Ra- 
dial velocity v r in Mollweide projection at 0.95Rq for case D5. This snapshot 
samples day 3880 (time tl) when the magnetic fields are strong, (b) Com- 
panion hydrodynamic case H5, whose stronger differential rotation shears 
out convective structures in the mid-latitudes, (c) Profile of mean angular 
velocity Q(r, 9) for case D5, with (d) radial cuts of Q at selected latitudes. 
(e) Meridional circulations for case D5, with magnitude and sense of circula- 
tion indicated by color (red counter-clockwise, blue clockwise) and stream- 
lines of mass flux overlaid. 

they are embedded. In the equatorial band, these structures 
move in a prograde fashion and at high latitudes in the retro- 
grade sense. In the polar regions, the radial velocity patterns 
have a somewhat cuspy appearance, with the strongest down- 
flows appearing to favor the westward and lower-latitude side 
of each convective cell. This may be a consequence of the 
strong retrograde differential rotation in those regions. 

The convective structures are quite similar to those real- 
ized in the hydrodynamic case H5 (Fig.[6j?), though there are 
some noticeable differences, particularly at the mid latitudes 
(around ±30°). In the hydrodynamic case there is little radial 
flow in these regions, as the strong differential rotation shears 
out the convective cells. This region is equatorward of the 
tangent cylinder, an imaginary boundary tangent to the base 
of the convection zone and aligned with the rotation axis. For 
rotating convective shells, it has generally been found that the 
dynamics are different inside and outside the tangent cylinder, 
due to differences in con nectivity and rotational constraint in 
these two regions (e.g.. lBusselll970l) . The tangent cylinder 
in our geometry intersects with the stellar surface at roughly 
±42° of latitude. In our compressible simulations, we gener- 
ally find that the convective patterns in the equatorial regions 



are bounded by a conic surface rather than the tangent cylin- 
der (Brown et al.l |20Q8|). In case H5 the strong differential 
rotation serves to disrupt the convection at the mid-latitudes. 
In contrast, in the dynamo case D5 the differential rotation 
is substantially weaker in both radial and latitudinal angular 
velocity contrasts (Table |2j. As is evident in Figure^, the 
convective cells fill in this region quite completely. 

The time-averaged angular velocity profile f2(r, 6) is shown 
for case D5 in Figures |6j:,d. The latitudinal angular velocity 
contrast AJ7i at and radial contrast Afl r in this case is remark- 
ably similar in amplitude to that realized in case D3 (Tabled, 
even though the basic rotation rate Slo is substantially faster. 
This is in marked contrast to our hydrodynamic companion 
cases where faster rotation leads to greater angular velocity 
contrasts. The accompanying meridional circulation patterns 
for case D5 (Fig. [6^) appear to have three major cells of cir- 
culation in each hemisphere. These flows are weaker than in 
case D3. They are very similar to the circulations found in 
case H5. 

5.2. Oscillations in Energies and Changes of Polarity 

The dynamo action realized by the convection in case D5 
exhibits significant changes in time. This time-varying be- 
havior is readily visible as oscillations of the volume-averaged 
kinetic and magnetic energy densities, as shown in Figure I7J? 
at a time long after the dynamo has saturated and reached 
equilibration. Here the kinetic energy of differential rotation 
(DRKE) undergoes factor of five changes on periods of 500- 
1000 days. As DRKE decreases the magnetic energies in- 
crease. Moving in concert are the mean toroidal (TME) and 
mean poloidal (PME) magnetic energies. The mean poloidal 
fields appear to lag slightly behind the mean toroidal fields as 
they both change in strength. The fluctuating magnetic en- 
ergies (FME) track the largest rises in the mean fields but 
decouple during many of the deepest dips. In contrast, the 
variations in convective kinetic energies (CKE) show little or- 
ganized behavior in time, and appear to change substantially 
only when the differential rotation is highly suppressed dur- 
ing the period from day 7500 to day 8300. The energy con- 
tained in the meridional circulations (MCKE) is weaker and 
not shown. Though it varies somewhat in time, there is not a 
clear relation to the changes in magnetic energies. 

Magnetic energies in case D5 can rise to be a substantial 
fraction of the kinetic energies. Averaged over the nearly 
16000 days shown here, the magnetic energies are about 17% 
of the kinetic energies. During individual oscillations the 
magnetic energies can range from a few percent of kinetic en- 
ergies to levels as high as 50%. The kinetic energy is largely 
in the fluctuating convection and differential rotation, with 
CKE fairly constant and ranging from 15-60% of the total 
kinetic energy as DRKE grows and subsides, itself contribut- 
ing between 40 to 85% of the kinetic energy. The magnetic 
energies are largely split between the mean toroidal fields and 
the fluctuating fields, with TME containing about 35% of the 
magnetic energy on average, FME containing about 61% and 
PME containing 4%. The roles of these energy reservoirs 
change somewhat through each oscillation. At any one time, 
between 10 and 60% of the magnetic energy is in TME while 
FME contains between 30 and 85% of the total. Meanwhile, 
PME can comprise as little as 1 % or as much as 10% of the to- 
tal. Generally, PME is about 12% of TME, but because PME 
lags the changes in TME slightly, there are periods of time 
when PME is almost 40% of TME. 

The global-scale magnetic fields can reverse their polarities 
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FIG. 7. — Complex time evolution in case D5 with flips in polarity of magnetic wreathes, (a) Volume-averaged energy densities for kinetic energy of convection 
(CKE), differential rotation (DRKE) and for magnetic energy in fluctuating fields (FME), in mean toroidal fields (TME) and in mean poloidal fields (PME) as 
labeled. Oscillations on roughly 500-1000 day periods are visible in the magnetic energies and in DRKE, though CKE stays nearly constant, (b) Mean toroidal 
field {Ba) averaged over entire northern and southern hemispheres (labeled) at mid-convection zone (0.85^0). Early in the simulation, opposite polarities 
dominate each hemisphere. Several reversals occur, along with several extreme excursions which do not flip the polarity of the global-scale field. During the 
interval from roughly day 7700 to 10200, the dynamo falls into peculiar single polarity states, with one polarity dominating both hemispheres. Bracketed interval 
from day 3500 to 5700 spans one full polarity reversal; (c) shows volume-averaged energy densities during this period, and (d) the mean toroidal field. Thick 
labeled tick marks above a,c indicate time samples used in later images. 



during some of the oscillations in magnetic energies. This 
is evident in Figure [7j? showing averages at mid-convection 
zone of the longitudinal magnetic field (B^) over the northern 
and southern hemispheres. Reversals in field polarity occur 
periodically, with typical time scales of roughly 1500 days. 
These reversals appear to happen shortly after peak magnetic 
energies are achieved, but do not occur every time magnetic 
energies undergo a full oscillation. Rather, it appears that sev- 
eral failed reversals occur where the magnetic energies drop 
and the average fields decline in strength, only to return with 
the same polarity a few hundred days later, for each successful 
polarity reversal. 

We focus in the following discussion on one such reversal, 
shown in closeup in Figures[7]:,fi and spanning the interval of 
time between days 3500 and 5700. Two reversals occur dur- 
ing this interval, with the global polarities flipping into a new 
state at roughly day 4100 and then changing back again at 
about day 5500. Detailed measurements of kinetic and mag- 
netic energies during this interval are shown in Table[3] 

5.3. Global-Scale Magnetic Reversals 

The nature of the global-scale magnetic fields during the re- 
versal spanning days 3500-5700 are presented in detail in Fig- 
ure [8] Several samples of longitudinal magnetic field B^ are 
shown at mid-convection zone spanning this time period. The 
timing of these samples is indicated in Figure [7] by numeric 
labels and likewise in Figure [8}? which shows azimuthally- 
averaged (B$) in a time-latitude map that spans the reversal. 

Before a reversal, the magnetic wreathes of case D5 are 
very similar in appearance to the wreathes realized in case D3. 
They are dominated by the azimuthally-averaged component 
of Z?0, while also showing small-scale variations where con- 
vective plumes distort the fields (Fig. |HJ»). At mid-convection 
zone, typical longitudinal field strengths are of order ±13 kG, 
while peak field strengths there can reach ±40 kG. Mean- 
while (B^) is fairly anti-symmetric between the northern and 
southern hemispheres (Fig. [8^). Shortly before a reversal, the 
magnetic wreathes strengthen in amplitude and become more 
anti-symmetric about the equator. 

They reach their peak values just before the polarity change 
at roughly day 4000 but then quickly begin to unravel, gain- 



ing significant structure on smaller scales (Fig. [5J;). At the 
same time, prominent magnetic structures detach from the 
higher-latitude edges and begin migrating toward the polar 
regions. Meanwhile, (B$) loses its anti-symmetry between 
the two hemispheres, with (B^) in one hemisphere typically 
remaining stronger and more concentrated than in the other 
(Fig.[8]'i). The stronger hemisphere (here the northern) retains 
its polarity for about 100 days as the fields in the other hemi- 
sphere (here southern) weaken and reverse in polarity. At this 
point the new wreathes of the next cycle, with opposite polar- 
ity, are already faintly visible at the equator. 

Within 100 days these new wreathes grow in strength and 
become comparable with the structures they replace, which 
are still visible at higher latitudes (Figs. [8}/,/). The mean 
(B<p) begins to contribute significantly to the overall structure 
of the new wreathes, and soon the polarity reversal is com- 
pleted. In the interval immediately after the reversal, small- 
scale fluctuations still contribute significantly to the overall 
structure of the wreathes, and B^ has complicated structure at 
mid-convection zone. At this time, the peak magnetic field 
strengths are somewhat lower, at about ±20 kG. As (B^) 
becomes stronger, the wreathes return to an anti-symmetric 
state, with similar amplitudes and structure in both the north- 
ern and southern hemispheres (Figs. [8^,y). They look much 
as they did before the reversal, though now with opposite po- 
larities. 

The wreathes from the previous cycle appear to move 
through the lower convection zone and toward higher lat- 
itudes. This can be seen variously in the time-latitude 
map at mid-convection zone (Fig. [SJj), in the Mollweide 
snapshots (Figs. [8]:-e) as well as in the samples of (B<p) 
(Figs. [HJi - j). This poleward migration is likely due to 
hoop stresses within the magnetic wreathes and an associ- 
ated poleward-slip instability (e.g.,ISpruit & van Ballegooiienl 
Il982t iMoreno-Insertis et alJll992h . Even at late times some 
signatures of the previous wreathes persist in the polar re- 
gions, and are still visible in Figures [8^,7 at day 4390. They 
are much weaker in amplitude than the wreathes at the equa- 
tor, but they persist until the wreathes from the next cycle 
move poleward and replace them. As they approach the po- 
lar regions, the old wreathes dissipate on both large and small 
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FIG. 8. — Evolution of longitudinal magnetic fields during a polarity re- 
versal, (a) Time-latitude plot of {B^) at mid-convection zone, (b—f) Snap- 
shots of in Mollweide projection at mid-convection zone (0.85K©) at 
times indicated by numbers in a. Between reversals the field is dominated by 
the mean component, but during reversals substantial fluctuations develop. 
(g— k) Accompanying samples of azimufhally-averaged (B^), showing struc- 
ture of mean fields with radius and latitude at same instants in time. 



scales, for the vortical polar convection shreds them apart and 
ohmic diffusion reconnects them with the relic wreathes of 
the previous cycle. 

Though reversals occur on average once every 1500 days, 
substantial variations can occur on shorter time scales. Here 
at mid-cycle the mean longitudinal field (B^) becomes quite 
weak as the wreathes become concentrated in smaller longi- 
tudinal intervals of the equatorial region (as in Figs. ^f,k at 
day 4780). At other times the mean longitudinal fields be- 
come quite asymmetric, with one hemisphere strong and one 
weak (i.e., during days 4900-5200) before regaining their anti- 
symmetric nature shortly prior to the next reversal. 

5.4. Temporal Changes in Differential Rotation 

The strong magnetic fields in case D5 suppress the global- 
scale flow of differential rotation. As the fields themselves 
vary in strength, the differential rotation responds in turn, be- 
coming stronger as the fields weaken and then diminishing as 
the fields are amplified. These cycles of faster and slower dif- 
ferential rotation are visible in the traces of DRKE shown in 
Figure [TJz. We revisit here the interval explored in closer de- 
tail both in Figure |7J: and in Figure [8] spanning days 3500 to 
5700 of the simulation and one full polarity reversal. 

The angular velocity ft at mid-convection zone is shown 
for this period as a time-latitude map in Figure^. Here again 
the timing marks tl-t5 refer to the snapshots of shown in 
Figures [SJ?- k. In the equatorial regions, the differential ro- 
tation remains fast and prograde, but with some modulation 
in time. Prominent structures of speedup are visible propa- 
gating toward the poles at the high latitudes. These structures 
are much more evident when we subtract the time-averaged 
profile of VL for this period at each latitude (Fig. [9j>). They 
appear as strong, tilted fast (red) structures extending from 
roughly ±30° latitude poleward. In the northern hemisphere, 
three such structures are launched over this interval. In con- 
trast, in the south only two such structures are evident. One 
is perhaps launched around day 4500 but does not survive or 
propagate. Comparing these features with the propagation of 
magnetic fields shown in Figure^ over the same interval, we 
find that velocity speedup features are well correlated with 
the poleward migration of mean longitudinal magnetic field. 
The velocity features bear some resemblance to the poleward 
branch of torsional oscillations observed in the solar convec- 
tion zone over the course of a solar magnetic activity cycle, 
though on a much shorter time scale here as befits the corre- 
spondingly shorter time between magnetic polarity reversals 
in these dynamo simulations. 

These velocity features propagate toward the poles rela- 
tively slowly. In a period of roughly 500 days they travel 
about 40° in latitude, or a distance of about 410 Mm. Their 
propagation velocity is about 0.8 Mm day -1 or about 9 ms -1 . 
This is considerably slower than the fluctuating latitudinal 
flows associated with the convection which at this depth have 
peak speeds of ±200 ms" 1 during this time period. The 
meridional circulations have amplitudes of about ±6 ms -1 
here but do not have a latitudinal structure at all similar to the 
pattern propagation. The propagation speed of the speedup 
patterns is closer to the Alfven velocity of the mean latitudi- 
nal magnetic fields, namely 

{B e ) 

V A,8 - n — - • V*i) 
V47Tp 

At mid-convection zone the mean density is about 
0.065 gem" 3 and {Be) in the poleward propagating plumes 
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FIG. 9. — Time-varying differential rotation in case D5. (a) Time-latitude 
map of angular velocity £2 at mid-convection zone (0.85Rq). There are sub- 
stantial temporal variations at both the equator and high latitudes, (b) These 
are accentuated by subtracting the time-averaged profile of Q(r. 9) at each lat- 
itude. Visible are poleward propagating speedup structures at high latitudes 
and more uniform modulations near the equator, (c) Corresponding varia- 
tions in Af2i at near the surface (upper curve, green) and at mid-convection 
zone (lower, red). 



ranges between roughly ±1.5 kG, yielding Alfven veloci- 
ties of about ±17 ms" 1 . This poleward migration may also 
be due to a poleward-slip instability arising from the strong 
toroidal fields. In this scenario, i f we neglect rotation and tur- 
bulent pumping in l atitude (e.g.. lMoreno-Insertis et alj|1992t 
iJouve & Brunl20 09). the propagation speed should be approx- 
imately the Alfven velocity of the mean toroidal magnetic 
fields, or about ±45 ms -1 based on a mid-convection zone 
Bcf, of approximately ±4 kG in the propagating features. If 
the poleward-slip instability is occurring, the velocity speedup 
features may result from conservation of angular momentum 
in the plasma that travels poleward with the wreathes. Rota- 
tion is likely to partially stabilize wreathes against poleward- 
slip dMoreno-In sertis et al1 l!992l) . and this may help explain 
their slower poleward propagation. The leaky topology of the 
wreathes will allow plasma to escape these structures, and 
this may modify the rate of their poleward propagation. The 
weak propagating features seen in case D3 (Fig. [5]l required 
nearly 700 days to propagate a similar distance in latitude. 
This difference may be due to the somewhat lower magnetic 
field strengths achieved in that dynamo. 

With the expanded sensitivity of Figure [9j?, we can see 
that the equatorial modulation appears as fast and slow pulses 
which span the latitude range of ±20°. These variations are 
fairly uniform across this equatorial region. The velocity vari- 
ations at the equator do not correspond with the equatorial 
propagating branch of torsional oscillations seen in the Sun 
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FIG. 10. — Extended history of temporally-varying differential rotation 
in Case D5. (a) Variations of Q(r,8) at mid-convection zone, (b) Tempo- 
ral variations are emphasized by subtracting the time-averaged profile of £7. 
Poleward propagating speedup structures are visible in each magnetic oscil- 
lation, (c) Variations in Af2] at near the surface (upper curve, green) and at 
mid-convection zone (lower, red). 



dThom pson et all [20031) . In the Sun, the equatorial branch 
may arise fr om enhanced cooling in the mag netically active 
regions (e.g.,[Spruit 2003; Rempel 200l l2007h . 

The temporal variations of the angular velocity contrast in 
latitude Af2i at is shown for this period in Figure [9J:. At mid- 
convection zone (sampled by red line) the variations in Af2i at 
are modest, varying by roughly 15%. Near the surface (green 
line) AQjat shows more substantial variations, with large con- 
trasts when the fields are strong in the magnetic cycle (prior 
to tl) and smaller contrasts when the fields are in the process 
of reversing (t2, t3). These near-surface values of Af2i at are 
reported in Table|2] averaged over this entire period (avg) and 
at points in time when the contrast is large (max, at day 3702) 
and small (min, at day 4060). 

6. SAMPLING MANY MAGNETIC CYCLES IN CASE D5 

The variations of angular velocities over considerably 
longer intervals of time for case D5 are shown in Fig- 
ures [l0h,fe. Here too we see the equatorial modulation over 
many magnetic cycles and the poleward propagating speedup 
bands. Asymmetries between the northern and southern hemi- 
sphere are evident at many times in different cycles. The lati- 
tudinal angular velocity contrasts shown in FigureflOb exhibit 
large variations. Successive magnetic cycles can have dis- 
tinctly different angular velocity contrasts, and there are addi- 
tional long-term modulations that span many magnetic cycles. 

A sampling of the associated magnetic field behavior is 
shown in the time-latitude maps of (B^) and (Bg) in Fig- 
ures Qjjj,/?. From day 1500 to 7300, four cycles occur in 
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FIG. 11. — Time-latitude plots for case D5. Shown at mid-convection zone 
are (a) mean longitudinal field (B^), and (b) mean colatitudinal field (Bg). 
Cycles of activity are visible, with fields changing polarity in the equatorial 
region. Also prominently visible are plumes of field reaching toward the p olar 
regions in a manner recalling Fig. [4] The time samples used in Figs.lliland ll2l 
are indicated. 

which wreathes of opposite polarity are achieved in each 
hemisphere. After this period, the dynamo explores unusual 
single-polarity states. Here either a single wreath is built (t7), 
or two wreathes of the same polarity (t6) occupy the two 
hemispheres. After day 10700 the dynamo emerges from this 
state and returns to building two wreathes of opposite polar- 
ity which flip in their sense an additional three times as the 
simulation continues. 

The poleward propagating magnetic features shown previ- 
ously in Figures [8}? are evident throughout this longer time 
sampling, now appearing as nearly vertical streaks in both 
(B^) and (Bg). These features continue to correlate well with 
the velocity speedup features evident in Figure fT0i>. 

6.1. Strange States and Wreathes of a Single Polarity 

These oscillating dynamos occasionally wander into dis- 
tinctly different states, and this occurs for case D5 around day 
7300. Instead of the two nearly anti-symmetric wreathes of 
opposite polarity above and below the equator, the dynamo 
enters a state where the polarity in each hemisphere is the 
same, as shown in Figures \Y2h.b at day 8903. Here two 
wreathes of same polarity occupy the two hemispheres and 
persist for an interval of more than 500 days. The positive- 
polarity Bff, reaches average amplitudes of 18 kG while the 
negative polarity structures have average amplitudes around 
3 kG. The azimuthally-averaged profiles of (B^) emphasize 
that these wreathes span the convection zone and have the 
same polarity everywhere. During this interval of time, the 
mean poloidal field has changed from an odd-parity state, with 
strong contributions from the odd-^ components, to an even- 
parity state where the even-£ components are more prominent. 

The dynamo can also achieve states where only a single 
wreath is built in the equatorial regions, as in Figures [12b. d at 
day 9590. Here a single strong wreath of positive polarity fills 
the northern hemisphere, with (Ba) reaching a peak amplitude 
of +18 kG. This unique structure persists for about 800 days 
before the dynamo flips polarity and builds a strong wreath 




FIG. 12. — Strange states in case D5. (a) Snapshot of B^ at mid-convection 
zone, showing two strong wreathes of the same polarity, (b) Instantaneous 
profile of (B^ ) at same time, (c) Snapshot of B^ at mid-convection zone at a 
time when a single wreath is formed, (d) Weaker negative polarity structures 
are visible in profile of (B^) at same instant, (e) Accompanying snapshot of 
IV at mid-convection zone, showing flows strongly affected by magnetism. 
(/) The instantaneous differential rotation, shown here as profile of f2(r, 0), is 
unaffected by the strong wreath. 

of negative polarity. The predecessor of this new wreath can 
be seen in profiles of (B$) where a much weaker structure of 
negative polarity is visible in the lower convection zone. 

The strong magnetic fields realized in the single wreath 
states react back on the convective flows. This is evident in the 
accompanying snapshot of radial velocities at mid-convection 
zone (Fig. fl2fe). In a narrow band spanning 0-20° latitude 
and coinciding with the strong tube, the upflows and down- 
flows have been virtually erased. Fluctuations in v<p and v$ 
are also very small in this region, and the flow is dominated by 
the streaming flows of differential rotation. Within the wreath 
the total magnetic energy (ME) at mid-convection zone is lo- 
cally about 10 to 100 times larger than the total kinetic energy 
(KE), while outside the wreath KE exceeds ME by factors of 
roughly 10 to 10 4 at this depth. We see similar restriction 
of the convective flows whenever the magnetic fields become 
this strong. 

The differential rotation itself (Fig. [l2]f ) is largely unaf- 
fected by the presence of the strong magnetic wreath. There 
is no clear signature of faster flow down the middle of the 
wreath. Likewise, there is no sign of the structure in profiles 
of the thermodynamic variables P, T, S, or p, with the mean 
profile instead dominated by latitudinal variations consistent 
with thermal wind balance. 

7. CREATING MAGNETIC WREATHES 

The magnetic wreathes formed in case D3 are dominated 
by strong mean longitudinal field components and show lit- 
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FIG. 13. — Generation of mean toroidal magnetic field in case D3. The view is from ±45° latitude to emphasize the equatorial regions, (a) Mean toroidal field 
(B^) with wreathes strongly evident, (b) Production by Pms serves to build (B^). This rate term generally matches the sense of (B^), thus being negative (blue 
in colorbar, with ranges indicated) in the core of the northern wreath and positive (red) in that of the southern wreath, (c) Destruction of mean toroidal field is 
achieved by the sum of the two fluctuating (turbulent) induction terms and the ohmic diffusion (Pps +Pfa +Pmd)- This sum clearly has opposite sense and similar 
magnitude to Pms- We break out these three destruction terms in the following panels, (d) Fluctuating (turbulent) shear Pps is strongest near the high-latitude 
side of each wreath, and (e) fluctuating (turbulent) advection Pfa is strongest in the cores of the wreathes. The sum of these terms (Pfs+Pfa) is responsible for 
about half the destructive balance, with the remainder coming from (/) the mean ohmic diffusion Pmd. Some differences arise in the boundary layers at top and 
bottom. 



tie variation in time. To understand the physical processes 
responsible for maintaining these magnetic wreathes, we ex- 
amine the terms arising in the time- and azimuth-averaged in- 
duction equation for case D3. 

7.1. Maintaining Wreathes of Toroidal Field 

We begin our analysis by exploring the maintenance of the 
mean toroidal field (B^). Here it is helpful to break the in- 
duction term from equation © into contributions from shear, 
advection and compression, namely 

V X(vxB) = 

(B-V)v-(v-'V)B- B(V_u) . (21) 

shear advection compression 

Details of this decomposition are given in the Appendix. 

The evolution of the mean longitudinal (toroidal) field (B^) 
is described symbolically in equation (1A8I >, with individual 
terms defined in equation JA91 >. When we analyze these terms 
in case D3, we find that (B^) is produced by the shear of dif- 
ferential rotation and is dissipated by a combination of tur- 
bulent induction and ohmic diffusion. This balance can be 
restated as 

« P MS + (Pfs + Pfa + Pvid) « , (22) 

at 

with Pms representing production by the mean shearing flow 
of differential rotation, Pps by fluctuating shear, Pfa by fluc- 
tuating advection, and Pmd by mean ohmic diffusion. Those 
terms are in turn 

Pms= ((B)-V)(v)\p (23) 
P FS = ((B'-V)v')\^ (24) 
P FA = -((«'-V)5%, (25) 
P M D = -VXT]VX(B)\,, (26) 



where brackets again indicate an azimuthal average and 
primes indicate fluctuating terms: v' = v - (v) . The detailed 
implementation of these terms is presented for our spherical 
geometry in equations (lAlOHAlST l. These terms are illustrated 
in Figure [13] for case D3, averaged over a 450 day interval 
from day 6450 to 6900. 

The structure of (B^) is shown in Figure [13b. The shear- 
ing flows of differential rotation Pms (Fig. \TSb ) act almost 
everywhere to reinforce the mean toroidal field. Thus the po- 
larity of this production term generally matches that of (B$). 
This production is balanced by destruction of mean field aris- 
ing from both turbulent induction and ohmic diffusion (sum 
shown in Fig. [T3b). The individual profiles of Pps, Pfa and 
Pmd are presented in turn in Figures[T3t/.e. f. The terms from 
turbulent induction (Pps and Pfa) contribute to roughly half of 
the total balance, with the remainder carried by ohmic diffu- 
sion of the mean fields (Pmd)- In the core of the wreathes, re- 
moval of mean toroidal field is largely accomplished by fluc- 
tuating advection Pfa (Fig. [13k ) and mean ohmic diffusion 
Pvid (Fig. [T3lf ), with the latter also important near the upper 
boundary. Turbulent shear becomes strongest near the bot- 
tom of the convection zone and in the regions near the high- 
latitude side of each wreath. Thus Pps (Fig.[T3t/) becomes the 
dominant member of the triad of terms seeking to diminish 
the mean toroidal field there. We find that the mean poloidal 
field is regenerated in roughly the same region. 

In the analysis presented in Figure Q~3] we have neglected 
the advection of (B$) by the meridional circulations (shown 
in the Appendix as PmaX which we find plays a very small 
role in the overall balance. We have also neglected the am- 
plification of (Be/,) by compressibility effects (the Appendix, 
Pmc and Pfc), though it does contribute slightly to reinforcing 
the underlying mean fields within the wreathes. 

To summarize, the mean toroidal fields are built through 
an il-effect, where production by the mean shearing flow of 
differential rotation (pyis) builds the underlying (B^). In the 
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FIG. 14. — Production of mean poloidal vector potential (A^) in case D3, 
with view restricted to ±45° latitude, (a) Mean poloidal vector potential 
(A$), with sense denoted by color (red, clockwise; blue, counter-clockwise). 
(b) The fluctuating (turbulent) emf £pi acts to build the vector potential. This 
term is strongest near the bottom of the convection zone and the poleward side 
of the wreathes, (c) Mean ohmic diffusion Emd acts everywhere in opposition 
to g pi. T he cores of the wreathes are positioned at roughly ±15° latitude 
(Fig.©). 



statistically steady state achieved, this production is balanced 
by a combination of turbulent induction (Pps +Pfa) an d ohmic 
diffusion of the mean fields (Pmd)- 

7.2. Maintaining the Poloidal Field 

The production of mean poloidal field is achieved through a 
slightly different balance, with turbulent induction producing 
poloidal field and ohmic diffusion acting to dissipate it. The 
mean flows play little role in the overall balance. This balance 
is clarified if we represent the mean poloidal field by its vector 
potential (A^), where 



(£? pol ) = (B r )f+ (Bg}9 = V X (A 



(27) 



as discussed in the Appendix. We recast the induction equa- 
tion (|4]l in terms of the poloidal vector potential by uncurling 
the equation once, obtaining 



dt 



(vxB)\,-r,Vx(B) 



(28) 



which is also equation dA29l ) in the Appendix. The first term 
is the electromotive force (emf) arising from the coupling of 
flows and magnetic fields, and the second term is the ohmic 
diffusion. These can be decomposed into contributions from 
mean and fluct uating components, as shown symbolically in 
equation ( 1A30I ). 

In case D3 we find that the mean poloidal vector potential 
(Afj,) is produced by the fluctuating (turbulent) emf and is dis- 
sipated by ohmic diffusion 



dt 



Epi + E 



MD ■ 



i0. 



(29) 



with £pi the emf arising from fluctuating flows and fluctuating 
fields, and contributing to the mean induction. The £md is the 
emf arising from mean ohmic diffusion. These terms are 



E m = (v'xB% = ( V %] 
Emd = — ?7 V X (B)\,. 



(v'B') 



(30) 
(31) 



The contribution arising from the omitted term Emi 
(see eq. IA31I ). related to the emf of mean flows and mean 
fields, is smaller than these first two by more than an order of 
magnitude. Additionally, Emi has a complicated spatial struc- 
ture which does not appear to act in a coherent fashion within 
the wreathes to either build or destroy mean poloidal field. 

The mean vector potential (A^) is shown in Figure [T4b . 
with poloidal field lines represented by the overlying con- 
tours. The mean radial magnetic field (B r ) is about ±1 kG 
in the cores of the wreathes, whereas the mean colatitudinal 
field (Bg) has an amplitude of roughly -2 kG (thus directed 
northward in both hemispheres), concentrated near the bottom 
of the convection zone. 

The production of (A^) by the fluctuating (turbulent) emf 
Epi is shown in Figure fl4"t >. Here too we average over the 
same 450 day interval. This term generally acts to reinforce 
the existing poloidal field, having the same sense as the un- 
derlying vector potential in most regions. It is strongest near 
the bottom of the convection zone and is concentrated at the 
poleward side of each wreath. This is similar, though not iden- 
tical, to the structure of destruction of mean toroidal field by 
fluctuating shear fps (Fig-HHf), It suggests that mean toroidal 
field is here being converted into mean poloidal field by the 
fluctuating flows. 

There are two terms that contribute to £pi, as shown in 
equation (f30b . Much of that fluctuating emf arises from corre- 
lations between fluctuating latitudinal flows and radial fields 
(— VgB' r ), which follows the structure of £pi (Fig. [14b) closely. 
The contribution from fluctuating radial flows and colatitudi- 
nal fields (v' r B' g ) is more complex in structure. Near ±20° 
latitude, this term reinforces (-v' e B' r ), but acts against it at 
higher latitudes and thus diminishes the overall amplitude of 
Em- The mean ohmic diffusion Emd (Fig. fT4b). almost en- 
tirely balances the production of (A$) by £pi- 

This shows that our mean poloidal magnetic field is main- 
tained by the fluctuating (turbulent) emf and is destroyed by 
ohmic diffusion. In mean-field dynamo theory, this is often 
parametrized by an "a-effect." Now we turn to interpretations 
within that framework. 

7.3. Exploring Mean-Field Interpretations 

Many mean-field theories assert that the production of mean 
poloidal field is likely to arise from the fluctuating emf. This 
process is often approximated with an a-effect, where it is 
proposed that the sense and amplitude of the emf scales with 
the mean toroidal field 



■/xB')=a(B) 



(32) 



where a can be either a simple scalar or may be related to 
the kinetic and magnetic (current) helicities. In isotropic (but 
not reflectionally symmetric), homogeneous, incompressible 
MHD turbulence 



at = -v' ■ (V X v') 



Olfn — 



4irp 



-B'-(VxB') 



(33) 
(34) 
(35) 



as 



discussed in iPouauet et all d 19761) and 
iBrandenburg & Subramanian (2005^ Here t is the life- 
time or correlation time of a typical turbulent eddy. In 
mean-field theory, these fluctuating helicities are typically 
not solved directly and are instead solved through auxiliary 
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FIG. 15. — Estimating the mean-field a-effect from case D3. (a) Fluc- 
tuating kinetic helicity a^. (b) Fluctuating magnetic (current) helicity a m . 
(c) Mean-field a, constructed by combining at and a m with a turbulent cor- 
relation time r. 

equations for the total magnetic helicity or are prescribed. 
Here we can directly measure our fluctuating helicities and 
examine whether they approximate our fluctuating emf. 

To assess the possible role of an a-effect in our simulation, 
we show in Figures [15b ,b the fluctuating kinetic and current 
helicities a* and a m realized in our case D3, averaged over the 
same 450 day analysis interval. To make an estimate of the a- 
effect, we approximate the correlation time r by defining 



where Hp is the local pressure scale height and v' is the local 
fluctuating rms velocity, which are functions of radius only. 
Estimated by this method, the turnover time r has a smooth 
radial profile and is roughly 10 days near the bottom of the 
convection zone, 3 days at mid-convection zone, and slightly 
less near the upper boundary. If we use the fast peak upflow 
or downflow velocities instead of the rms velocities, our es- 
timate of t is about a factor of 4 smaller. Our mean-field a 
(eq. [33j is shown in FigurefTBl:. In the upper convection zone, 
this is dominated by the fluctuating kinetic helicity while the 
fluctuating magnetic (current) helicity becomes important at 
depth. 

We form a mean-field emf (right-hand side of eq. [32l > by 
multiplying our derived a (Fig. [TBI:) with our (B^) (Fig.QjJj), 
and show this in Figure [16b . The turbulent emf Ejh, which 
is the left-hand side of equation d32l >. can be measured in 
our simulations and is shown again in Figure [T6i>. Although 
there is some correspondence in the two patterns, there are 
significant differences. In particular, the mean-field emf has 
peak amplitudes in the cores of the wreathes (at ±15° lati- 
tude) and is negative there. In contrast, the actual fluctuat- 
ing emf given by £fi is positive and has its highest ampli- 
tude at the poleward side of the wreathes. Thus the mean- 
field emf predicts an incorrect balance in the generation terms 
and would yield a distinctly different mean poloidal mag- 
netic field. To assess whether better agreement may be 
achieved with a latitude-averaged emf, we average the mean- 
field emf and £pi separately over the northern and southern 
hemispheres and plot these quantities in Figure [T6l:. Though 
both have a similar positive sense near the base of the convec- 
tion zone, the hemisphere-averaged £pi becomes small above 




FIG. 16. — Comparison of emfs in case D3. (a) Profile of proposed mean- 
field emf given by aft^). (b) Actual turbulent emf Ef\ measured in the 
dynamo, (c) Variation of hemisphere-averaged emfs with fractional radius. 
The mean-field approximated emf is shown in blue, and £pi in red. The 
average over the northern hemisphere is shown solid, the southern is dashed. 



O.87?0 whereas the averaged mean-field emf a(B ( f > ) is large 
and negative there. Thus even the averaged emfs are not in 
accord. 

In summary, it is evident that a simple scalar a-effect will 
predict the wrong sign for the fluctuating emf in the two hemi- 
spheres, as (B$) is anti-symmetric across the equator while 
(A^) is symmetric. An a-effect based on the kinetic helicity 
and magnetic helicity may capture some sense of the fluctu- 
ating emf, as those quantities are themselves anti-symmetric 
across the equator. Yet Figure [16] suggests that there are 
significant discrepancies between this particular approxima- 
tion and our turbulent emf. In particular, this mean-field a- 
effect misses the offset between the generation regions for 
mean toroidal and mean poloidal field. This offset in lati- 
tude of the generation regions may be important for avoid- 
ing the a-quenching problems encountered in many mean- 
field theories. A more complex mean-field model, which 
takes spatial gradients of (B^) into account, may do better - 
In particular, the O X J-effect ( e.g.,| Moffatt & Proctori ll982t 
Rogachevskii & Kleeorinll2003l) may be at work in these sys- 
tems, and preliminary explorations indicate that this term 
matches the spatial structure of our £pi better than the above 
a-effect. A tensor representation of the a-effect may do much 
better, and test-field te chniques could be em ployed to mea- 
sure this q uantity (e.g.,[Sc hrinner e t al.ll2005l and recently re- 
viewed in Brandenbur gl2009l) . 

8. CONCLUSIONS 

The ability for a dynamo to build wreathes of strong mag- 
netic fields in the bulk of the convection zone has largely 
been a surprise, for it had generally been supposed that tur- 
bulent convection would disrupt such magnetic structures. 
To avoid these difficulties, many solar and stellar dynamo 
theories shift the burden of magnetic storage, amplification 
and organization to a tachocline of shear and penetration at 
the base of the convection zone where motions are more 
quiescent. In contrast, our simulations of rapidly rotating 
stars are able to achieve sustained global-scale dynamo action 
within the convection zone itself, with the magnetic structures 
both being built and able to survive while embedded deep 
within the turbulence. These dynamos are able to circum- 
vent the Parker instability by means of turbulent Reynolds and 
Maxwell stresses that contribute to the mechanical force bal- 
ance and prevent the wreathes from buoyantly escaping the 
convection zone. This striking behavior may be enabled by 
the stars rotating three to five times faster than the current Sun, 
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which yields a strong differential rotation that is a key element 
in the dynamo behavior. It is quite interesting that in our dy- 
namo cases the angular velocity contrast in latitude and radius 
is almost constant at differing rotation rates, whereas our hy- 
drodynamic cases tend to have increasing AO with more rapid 
rotation ilo- 

We have achieved some dynamo states that are persistent 
and others that flip the sense of their magnetic fields. In our 
case D3 the global-scale fields have small vacillations in their 
amplitudes, but the magnetic wreathes retain their identities 
for many thousands of days. This represents hundreds of rota- 
tion periods and several magnetic diffusion times, indicating 
that the dynamo has achieved a persistent equilibrium. In- 
creasing the rotation rate yields more complex time depen- 
dence. In our case D5 the oscillations can become large, and 
may result in the global-scale fields repeatedly flipping their 
polarity. At times this dynamo appears to be cyclic but in 
other intervals it behaves more chaotically. 

In our persistent case D3 we are able to analyze the gener- 
ation and transport of mean magnetic field. We find that our 
dynamo action is of an a - nature, with the mean toroidal 
fields being generated by an £!-effect from the mean shearing 
flow of differential rotation. This generation is balanced by a 
combination of turbulent induction and ohmic diffusion. The 
mean poloidal fields are generated by an a-effect arising from 
couplings between the fluctuating flows and fluctuating fields, 
with this production largely balanced by the ohmic diffusion. 
This is unlike the toroidal balance, for here the mean flows 
play almost no role and the turbulent correlations are con- 
structive rather than destructive. In assessing what a mean- 
field model might predict for the magnetic structures realized 
in case D3, we find that the isotropic, homogeneous a-effect 
based on kinetic and magnetic (current) helicities fails to cap- 
ture the sense of our turbulent emf. In general, our Epi is 
poorly represented by an a(B^) that is so determined. 

The transition to richer time dependence with increasing ro- 
tation rate in case D5 appears to arise from subtle changes and 
phasing relationships between the toroidal and poloidal mag- 
netic production terms. These are difficult to assess, since 
the production terms themselves are complex in space and 
now vary in time as well. It is evident that the mean poloidal 
fields lag the changes in the mean toroidal fields, and there is 
clear migration of magnetic field to the higher latitudes dur- 
ing a reversal. This latitudinal migration may result from a 
polarward-slip instability, triggered by the stronger magnetic 
fields that are generated by this dynamo, and this migration 
may lead to the oscillations in field strength and polarity. The 
analysis of magnetic field production that we have carried out 
for case D3 required significant averaging of the turbulent 
correlations, accomplished here by time averages, to obtain 
a coherent picture of the balances achieved. This has not yet 
been tractable for the oscillating solutions, as the generation 
terms change on shorter time scales than needed to obtain sta- 
ble averages of the turbulent processes. It appears that the 
imbalances in the production and destruction of mean mag- 
netic fields during cyclic behavior are modest, and currently 
we cannot pinpoint just which terms out of the large medley 
serve to drive the oscillations. 

These dynamo oscillations are not special to case D5. In- 
deed, we have explored a broader class of oscillating dynamo 
solutions, which will be reported in a subsequent paper. Some 
of these solutions are realized by taking our more slowly ro- 
tating case D3 to higher levels of turbulence by reducing the 



eddy diffusivities, while others are achieved at even higher 
rotation rates. We find such global-scale oscillations and po- 
larity reversals fascinating, since these appear to be the first 
self-consistent 3-D stellar dynamo simulations which achieve 
such temporally organized behavior in the bulk of the convec- 
tion zone. 

Accompanying the oscillations in global-scale magnetic 
field are changes in the differential rotation, and these sig- 
natures are the most likely to be found through stellar ob- 
servations. The angular velocity contrast Af2i at between the 
equator and high latitudes can vary by 20% over periods of 
hundreds or thousands of days. The patterns of speedup 
which we find propagating toward the polar regions of these 
dynamos may have correspondence with the polar branch 
of the torsional oscillations which are observed in the Sun. 
The speedup features may arise from conservation of angu- 
lar momentum in fluid which is partially trapped within the 
wreathes as they slip toward the poles. The leaky nature of 
our wreathes will modify this somewhat and that may ex- 
plain why the wreathes are not stabilized by rotation (e.g., 
iMoreno-Insertis et al.ll 19921) . 

The realization of global-scale magnetic structures in our 
simulations, and their great strength relative to the fluctuating 
fields, may in part be a consequence of the relatively modest 
degree of turbulence attained here. Whether such structures 
can be generated and sustained amidst the far more complex 
flows in actual stellar interiors is not yet clear. If such struc- 
tures are indeed realized in stars, they may or may not sur- 
vive to print through the highly turbulent convection occur- 
ring just below the stellar photosphere. If they do appear at 
the surface, some global-scale magnetic features may propa- 
gate toward the poles along with the bands of angular velocity 
speedup. There are some indications in stellar observations 
that global-scale toroidal magn etic fields may indeed become 
strong in rapidly rotating stars dDonati et al.ll2006t iPetit et all 
2008), though small-scale fields may still account for much o f 
the magnetic energy near the surface dReiners & Basrill2009h . 
The global-scale poloidal fields may be more successful in 
surviving the passage through the turbulent surface convec- 
tion. If they do, the stellar magnetic field will likely have 
significant non-dipole components. Thus the mean poloidal 
fields observed at the surface may give clues to the presence 
of large wreathes of magnetism that occupy the bulk of the 
convection zone. 
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APPENDIX 

PRODUCTION, DESTRUCTION AND TRANSPORT OF MAGNETIC FIELD 

We derive diagnostic tools to evaluate the generation and transport of magnetic field in a magnetized and rotating turbulent 
convection zone. This derivation is in spherical coordinates, and is under the anelastic approximation. 

Induction Equation 

In the induction equation the first term on the right hand side represents production of magnetic field while the second 
term represents its diffusion. We first rewrite the production term to make the contributions of shear, advection and compressible 
effects more explicit as 

V X (v X B) = (B ■ V)v-(v ■ V)B-B(V-v). (Al) 

Under the anelastic approximation the divergence of v can be expressed in terms of the logarithmic derivative of the mean density 
because 



and therefore 



The induction equation thus becomes 



V- (pv) = = p( V- v) + (v ■ V)p, 



d 

V-y = -v r — In p. (A2) 
Or 



— = (B ■ V)v- (v ■ V)B + v r B— \np - V X (77V X B) (A3) 

shearing advection v y . diffusion 
b compression 

As labeled, the first term represents shearing of B, the second term advection of B, the third one compressible amplification of 
B, and the last term ohmic diffusion. 

Production of Axisymmetric Magnetic Field 

To identify the processes contributing to the production of mean (axisymmetric) field, we separate our velocities and magnetic 
fields into mean and fluctuating components v = (v) + v' and B = (B) + B' whe re an gle brackets denote an average in longitude. 
Thus (v 1 ) = (B 1 ) = by definition. Expanding the production term of equation (I A3t we obtain the mean shearing term 



the mean advection term 



and the mean compressibility term 



((B ■ V)v) = ((B)- V) (v) + ((£?' • V)v') , (A4) 
- ((» • V)S) =-((«)• V) (B) - (W ■ V)B') , (A5) 



(v r B— lnp) = ({v r ){B) + {v' r B')) —\np. (A6) 
or v 'Or 



In a similar fashion, the mean diffusion term becomes 

-(Vx(i)VxB))=-Vx(i)Vx(B». (A7) 
The axisymmetric component of the induction equation is written symbolically as: 

d(B) 



Of 



= Pms + Pfs + Pma + Pfa + Pmc + Pfc + Pmd ( A8) 



With Pms representing production of field by mean shear, Pps production by fluctuating shear, Pma advection by mean flows, Pfa 
advection by fluctuating flows, Pmc amplification arising from the compressibility of mean flows, Ppc amplification arising from 
fluctuating compressible motions, and Pmd ohmic diffusion of the mean fields. In turn, these terms are 

Pms = ((B)-V)(v), P fs =((B'-V)v'}, P ma =-((v)-V)(B), P fa =-((«'• V)B'), 

d d 

Pmc = ((v,.)(B))— lnp, P F c = ((v;B'))^:lnp,and P MD =- V X (77V X (B)). (A9) 
We now expand each of these terms into their full representation in spherical coordinates. 
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Production of Mean Longitudinal Field 



dt 

Pms - 
Pfs - 
Pma z 
Pfa~- 



-Pms + Pfs + Pma + Pfa + Pmc + ^fc + ^md 

(B^(v r )+cot9(B^)(vg) 



(v*> + - 



r dr r dO rsin 



V rl dr r 89 



(B 6 )- 



{v^)(B r ) + cote{ V<k ){Be) 



' dr r 86 rsin 9 defy 
d 



{v'B' r ) + C oi9{v'B' e ) 



Pmc= ((v r )(B^,)) — lnp 



P FC = {{v' r B' 4> ))-lnp 



Pmd= VV 2 (B, 



r?(B ) rfr? /15(r(fi )) 



r 2 sin 2 6> dr \r dr 
Production of Mean Latitudinal Field 



(A10) 
(All) 
(A12) 
(A13) 
(A14) 
(A15) 



d(B e ) 



dt 

Pms = 



= Pms + Pfs + Pma + Pfa + Pmc + Pfc + Pmd 

(B )(v r )-cot9(B (P }(v^ 



{Br) dr + —~de 



(ve} + - 



Pfs = 
Pma = ~ 
Pfa = ~ 



Bl 0B A 8 B^ d 



r dr r d9 rsm9 dd> 



vo ) + 



{B' g v' r )-cot9{B'^) 



dr r d9 

, d v'a d 
v' — + + ■ 



(Be)~ 
v's d 



{v e ){B r )-cote{v 4 ,){B 4 



* dr r d9 rsm9 dip 
d 



B' n 



(v e B r )-cot6(v'B>) 



Pmc= {{v r ){B 6 ))—\np 



^FC '■ 



(W))|lnp 



Pmd= V^ 2 (Bg) + 



2 v d(B r ) v (B e ) dr, ( \ d(r{B e )) \d{B T ) 



r 2 d9 r 2 sm 2 9 dr \r dr r d9 
Production of Mean Radial Field 



(A16) 
(A17) 
(A18) 
(A19) 
(A20) 
(A21) 



d(B r ) 



dt 

Pms = 



- Pms + Pfs + Pma + Pfa + Pmc + Pfc + Pmd 

(Bg)(vg) + (B 4> }(V 4> } 



(Vr>- 



Pfs = 
Pma = ~ 
Pfa=~ 



' dr r d9 rsin 



{Vr) o + Md 



(B r } + 



(v e }(Bg} + (v 4> }(B 4> ) 



dr r d9 
,d v' e d ^ d 



" dr r d9 rsm9 defy 



b:)+ 



d 



Pmc= ((v r )(B r )) —Inp 



P FC =((vX))-lnp 



Pmd= 7?V z (B r )-27?^ 



B r ) 2r]d{B e ) 2r]cot9(Bg) 



(A22) 
(A23) 
(A24) 
(A25) 
(A26) 
(A27) 
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Maintaining the Poloidal Vector Potential 

The balances achieved in maintaining the mean poloidal magnetic field are somewhat clearer if we consider its vector potential 
rather than the fields themselves. The mean poloidal field (B po \) has a corresponding vector potential (A^), where 

(B poi ) =(B r )r+(B e )e = Vx(A\ 4> ) 

= ThTk( A <^ e )?-7T r ^)Q (A28) 

= V X /a </> 



The other components of the poloidal vector potential disappear, as terms involving d/d(f) vanish in the azimuthally-averaged 
equations. Likewise, the 0-component of the possible gauge term VA is zero by virtue of axisymmetry. We recast the induction 
equation (eq.|4]l in terms of the poloidal vector potential by uncurling the equation once and obtain 

^M=vxB\ <f> -r ] VxB\ <l> . (A29) 

This can then be decomposed into mean and fluctuating contributions, and represented symbolically as 

d(A ) 



dt 



= E M i + E F i + E M d, (A30) 



with Emi representing the electromotive forces (emf) arising from mean flows and mean fields, and related to their mean induction. 
Likewise, £pi is the emf from fluctuating flows and fields and Emd is the emf arising from mean diffusion. These are in turn 

E Ml = (v)x{B)\ <f> = {v r ){B e )-{v e ){B r ), (A31) 
E m = (v'xB')\ (j> =( V ' r B' e )-(v'eK), (A32) 
E MD = - V V X (B)\,=-V- r (A33) 

Fluctuating (Non-Axisymmetric) Component of the Induction Equation 

Left out of this analysis is the fluctuating component of the induction equation. This can be derived by subtracting the mean 
induction equation (IA8b from the full induction equation, yielding the following equation for the fluctuating fields 

BB' 

— = ((B)-V)v' + (B'-V){v)+e 
-i{v}-V)B'-(v'-V){B)-T 

+((v r )B' + v' r (B))^-lnp+g 
or 

-V X (77V X (B')) (A34) 

where the quantities £ = (B' ■ V)v' - ((B' ■ V)v'), T = (V • V)B'-((v' ■ V)B'), and Q = (v' r B' - {v' r B'))-§-\np, represent the 
difference between mixed stresses from which we subtract their axisymmetric mean. In the standard mean-field derivation, these 
quantities are siblings of the G -current involv ing the mean electromotive force (v X B) and its 3-D equivalent v X B (i.e., the so 
called "pain in the neck" term, lMoffattlll978"b . 
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